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ABSTRACT 

We develop a new approach for studying flux anomalies in quadruply-imaged fold lens systems. 
We show that in the absence of substructure, microlensing, or differential absorption, the expected 
flux ratios of a fold pair can be tightly constrained using only geometric arguments. We apply this 
technique to 11 known quadruple lens systems in the radio and infrared, and compare our estimates to 
the Monte Carlo based results of Keeton, Gaudi, and Petters (2005). We show that a robust estimate 
for a flux ratio from a smoothly varying potential can be found, and at long wavelengths those lenses 
deviating from from this ratio almost certainly contain significant substructure. 

Subject headings: gravitational lensing - galaxies: structure - galaxies: fundamental parameters 
(masses) 
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1. INTRODUCTION 

To date, at least forty-eight multiply lensed quasars 
with three or more images have been discovered 
(e.g. CASTLES 1 , IHewitt et al.lU^9l llnada et all [20051: 
iMore et al]|2009ft . The great advantage to studying mul- 
tiply imaged systems is that symmetries of the lensing 
galaxy can b e understood witho ut a detailed mass re- 
construction (jPetters et alj|200l . 

Of particular int erest are the so-called "fold" lenses 
(jKeeton et al.ll2005l hereafter KGP) in which two images 
lie on opposite sides of the tangential critical curve in the 
image plane, while in the source plane, the source lies 
near an edge of a tangential caustic. We illustrate the 
geometry of a fold lens in Fig. [1] 

As a source gets arbitrarily close to the caustic it 
can be shown from purely analytic arguments that mag- 
nificati on of the two images should be equal and op- 
positc (Blandfo rd fc Naravanl [l98l IPetters et all [2001: 
ISchneider et all 11992). Thus, we would naively expect 
the fold relation: 



Rfoid 



fA - Jb 
Sa + Sb 



(1) 



to be zero, where /a,b are the fluxes of the individual 
images in the same band. Our convention is that images 

with a negative parity still hav e a positive flux. 

Observationally, h owever (e.g. lPoolev et al.l2006l 120071 
IKeeton et al.l I2006D there is often a significant flux 
anomaly between two images, there has been signifi- 
cant discussion on the nature of the fold flux anoma- 
lies (IMao fe Schneiderl 119981 IDalal fe Kochanekl [200l 
ICongdon fc Keetonll2005rr For optical lenses, two of the 
most common explanatio ns include: microlensing from 
stars in the lensing galaxy (jK oopma ns fe de BruviJl2000 ; 



Metcalf & Madau 2001; Schechter & Wambsganss 2002 
Keetonl 120031: iChartas et al.l 12004 iMorgan et al.l 12006 



Anguita et al.l I2008D a nd differential reddening by dust 
(|Lawrence et al.l I1995T ). These causes of flux anomalies 
are expected to be highly wavelength-dependent, how- 
ever. For example, differential absorption will strongly 
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affect optical photometry, but will have limited or no ef- 
fect in the infrared (IR) or radio. Though microlensing 
applies achromatically to the lensed image, it is only a 
significant effect when the Einstein radius of the lens- 
ing star is similar to or larger than the angular size of an 
emission region in a particular waveband. Typically, Ein- 
stein radii of individual stars at cosmological distances 
(i.e., in the lensing galaxy) will be on the order of mi- 
croarcseconds which corresponds to a typical angular size 
of quasar optical emission regions. However, radio emis- 
sion, especially from radio lobes, is generally much more 
extended. As a result, microlensing is less likely to cause 
significant flux anomalies in the radio. 

While the current work focuses exclusively on fold 
lenses, f uture geometric analysis may shed ligh t on "cusp 
lenses" (jKeeton et al.ll2003b IPetters et al.ll2001h . in which 
a source image near the corner of a caustic produces 3 
clustered images in the foreground plane. As with folds, 
cusps are naively expected to obey a simple flux ratio re- 
lation in which the brightest image precisely equals the 
sum of the fluxes of the dimmer two. As with folds, 
significant deviations from this expectation have been 
observed. 

Even if we confine the discussion to observations 
at long wavelength, fold flux anomalies don't disap- 
pear. This is likely attributable t o small-scale varia- 
tions in the lens galaxy po tential (|Congdon fc Keetonl 
l2005HMao fc Schneiderlll998| ). To that end, most workers 
have focused on generating semi-analytic models which 
closely fit the observed image positions and fluxes. They 
propose adding explicit substructure to their models as 
a kludge to correct the fluxes. There has been enor- 
mous eff ort expended trying to model galaxy le nses ex- 
plicit ly (jMunoz et al.l [2001) using lensmodel (Keeton 
12001ft and other software. However, for many of t hese 
lenses (|Evans fc Wittl l200l IMao fc Schn"eide^ll998l ) no 
simple analytic model will suffice. Some res earchers 
(|Evans fc Wittl l200l ICongdon fc Keetonl I2005T) decom- 
pose lensing potentials into an orthonormal basis sets or 
multipole moments and show that any set of observa- 
tional constraints may be fit with a complex enough ex- 
pansion of the potential. It is noteworthy, however, that 
all of these authors stress that not all of these models 
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provide viable explanations of the flux anomalies. 

In reality there is a great deal of information about 
the local lensing field which can be gleaned geometri- 
cally. That is, using only the positions of the observed 
lensed images, we will derive a semi-analytic approach 
to smooth lens flux anomalies. A particular system can 
then be analyzed without recourse to complicated mod- 
els. 

In this paper, we will develop a semi-analytic geomet- 
ric approach to understanding fold lenses. This approach 
is intended as a first step toward a more general theory 
involving measurement of galaxy lens substructure. Our 
approach is as follows. In ^ we derive a geometry-based 
expansion of a smooth potential, and discuss how observ- 
ables can be used to uniquely predict a flux anomaly. In 
$3]we test our semi-analytic results with simulated galaxy 
lenses, and show that our approach allows us to identify 
lenses with significant substructure. In $4l we describe 
an observational set of fold lenses in the radio and IR 
regime, and in $5j we apply our sample. Finally, in S}6] 
we discuss future prospects. 

2. THEORY 
2.1. Notation 

It is useful to give a bit of background regarding the 
notation we'll be using throughout. As is the normal 
practice, we will use the angular vector, (3 to describe 
(non-observable) positions in the background plane in the 
absence of lensing. Likewise, we use 6 for the observed 
position(s) of the lensed image. They are related via: 



p i = O i -a i {0) 



(2) 



where i = 1,2 used for the principle directions, and the 
displacement vector is defined as: 



at = VV* 



(3) 



The subscripted index represents a single angular deriva- 
tive in the 9i direction. Since we will be performing mul- 
tiple derivatives, subsequent derivatives of the potential 
will omit the comma in order to reduce clutter. 

The potential is generated via the two-dimensional 
Poisson equation: 

V 2 V = 2k (4) 

where /t, the convergence, is the dimensionless surface 
density of the galaxy. 

Likewise, the Jacobian of the source position generates 
two shear terms: 



7i : 



V>U - V>2 



where 



72 = V'12 



2 2,2 

7 = 7i + 7 2 



(5) 
(6) 
(7) 

(8) 



This yields an inverse magnification field of: 

M - 1 = (l- K ) 2 -7 2 (9) 
where along "critical curves" of the lens, = 0. 



2.2. Fold Flux Anomalies 

KGP derived an analytic expression for the expected 
fold relation arising from a smooth potential. Their ex- 
pression involved the Taylor expansion of the potential 
around a point on the critical curve. In order to esti- 
mate the fold relation, Rfold, in the simplest form, they 
performed a rotation around the center of the lens such 
that the fold images are oriented vertically from one an- 
other, with the positive par ity image ("A") below the 
negative parity one ("B"). In lKochanek et all (|2004l ). the 
interested reader can find some very helpful figures to il- 
lustrate image parity. In Fig. [2] we show our rotated 
coordinate system (along with a few other angles to be 
used later). 

By definition, along the critical curve, equation ([9]) 
equals zero. The rotated coordinates were specifically 
selected such that V'12 = 0, and V22 = 1- This rotation 
can be performed without a loss of generality. Likewise, 
for our derivation, we further allow a reflection such that 
the rotated images appear in the positive-x half of the 
plane. 

We may parameterize the fold ratio with a simple form: 

Rfold — Af iddi (10) 

where throughout, we will refer to Af Q id as the "anomaly 
parameter", and Rfold as the "fold relation." The 
anomaly parameter is introduced because it can be shown 
to be a constant in the limit of small separations. KGP 
showed that the anomaly parameter may be expressed 
as: 



A 



fold 



3V>?22 



3^112^222 + ^2222(1 - V^l) 



6^222(1 - -011) 



(11) 



and d\ is the distance between the two images in the 
observed plane. We note that this expression contains 
one term (V'ii) related to the local shear field, many (the 
third derivatives) related to the local flexio n field (e.g., 
IGoldberg fc Bacodl200l iBacon et afll2006| ), and 1 term 
related to the fourth derivative, which is normally not 
considered at all in lensing analysis. 

2.3. Simplified Analytic Models 

KGP model a number of observed multiple lensed sys- 
tems using a wide range of elliptical isothermal models, 
and show that there is a small reasonable range of flux 
anomalies that might be expected for known fold sys- 
tems. However, for smooth potential models, simulations 
aren't necessary. For this paper, we will define a "smooth 
potential" very strictly, but will subsequently show that 
even non-smooth potentials can be adequately fit. For 
our derivation, we assume: 

1. That the potential must be expressible as a circu- 
larly symmetric potential, plus an external shear: 

V/ tot ) (0) = V(6>) + 7 e 2 cos(2z/) (12) 

where j e is the magnitude of the external shear, 
and v is the angle between the induced shear and 
the radial vector. In our rotated coordinate frame, 
this is essentially the angle that the best-fit critical 
curve ellipse makes with the horizontal. 

2. That the shape of the critical curve is largely inde- 
pendent of a detailed model of the potential. This 
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is due in part to the fact that the critical curve 
must thread between the 4 images in a quad sys- 
tem. This is one of the main reasons that the image 
positions of quad systems may be fit very generi- 
cally while fluxes are more complicated to fit. 

In other words, we will generate a simple lens 
model using an isothermal sphere plus external 
shear based on positions only. We argue that the 
shape of the critical curve will not vary significantly 
from other models, an assumption we will test in 



3. All of our expansion will be around an as yet un- 
known point, V, which lies along the critical curve. 
This is the same point used in equation (jlip to de- 
termine the A f id parameter (and related to that, 
the local derivatives of the potential). The vector 
between the center of the lens and V makes an an- 
gle, r\ with respect to the horizontal (as illustrated 
in FigEJ) . This angle is assumed to be small. For 
the sample discussed in fjH the maximum 77 is 26 
degrees with an average over the sample of only 12 
degrees. 

With full generality, we can apply local geometric prop- 
erties of the critical curve to simplify equation (|TT]) . 
First, we note that at V: 

v^- 1 ) = 

[-011l(l ~ 022) - "0122(1 - 0n) ~ 2^12^112] "i 
+ [-'0112(1 - ?/>22) - "0222(1 - 1p n ) - 2^120122] j 
= -(1 -0il)(0122« + 0222j) (13) 

where the simplification can be found by applying the 
constraints on the second derivatives from the choice of 
rotation. 

Since the gradient is perpendicular to the critical curve, 
if the curve makes an angle, 0, with vertical (as shown 
in FigEJ then: 

0222 (14) 



tan( 



(16) 



tan m — 
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again, with complete generality. 

We now Taylor expand the circular component of the 
potential around V, such that: 

2 b 24 

(15) 

where V is a distance 9 from the center of the lens and 
primes represent radial derivatives of the circular compo- 
nent of the potential. We can then expand equation fj 15|) 
explicitly, such that X = 9o cos 77 and 8 y — 9q sin 77. 

For an assumed value of 77, we can thus compute all 
higher derivatives of the potential. While the Taylor ex- 
pansion does not include the external shear component, 
it should be noted that the only term in equation (TTT1) 
which will be affected by the external shear is 0n. All 
3rd derivatives and higher will include only the circular 
component of the potential. 

Thus, we may expand the ratio as a series in sin 77 



0-222 „ . 1 

- = 3 sin 77+ - 

0122 2 V 0' 



150"6» O 



10" 



150' \ 

sin 77 + 



For an isothermal circular component, the quantity in 
the parentheses reduces to 15, while for a point mass, it 
becomes 4. In any event, under under the assumption 
that 77 is small, only the first term matters, and we get a 
relationship which is largely independent of radial profile. 
For small 77: 

~ -377 (17) 

Thus, for a fiducial model of the critical curve near 
the fold, a unique position, V, can be found. As a first 
step in the procedure, consider the critical curve near the 
midpoint of the two fold images. About the midpoint, 
measure the local curvature. Using standard trigonom- 
etry, this arc uniquely defines a circle in the lens plane 
with a center at coordinates (x c , y c ), with radius of cur- 
vature, r, from which some tedious algebra yields: 



sm 77 



Vc 



r sin < 



(18) 



where 



C — \fr 2 + x 2 . + y 2 + 2x c r cos < 



(19) 

Combined with equation (|17j) . 77 can be solved iteratively 
very quickly. 

Objections might be raised that application of equa- 
tion (fT5)) requires that we have a global model of the 
potential. This is true, but as we will show in §[3j a wide 
range of models will leave the local shape of the critical 
curve largely unchanged. Moreover, for most lenses, a 
good estimate of 77 can be found by simply taking the 
midpoint of the two fold images. 

2.4. Semi-Analytic Fold Ratio Estimates 

Once we have an estimate of 77, it is a straightforward 
matter to estimate the flux anomaly parameter, Af a id 
in terms of the radial derivatives of the potential field 
using equation (|f 5[) . Each of the terms in equation (jlip 
may then be expanded explicitly. If we define a smooth 
field as a circular profile plus an external shear, only 2nd 
derivatives contain any indication of the non-circularity. 

In particular, note that to first order in 77: 



0: 
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^' /n ^ 1 

7 cos(2i^) = 1 



(20) 



where the last equality is guaranteed by the choice of 
coordinates. Thus, we have: 



<tf = 6>o [1 + 7 e cos(2i/)] 



(21) 



It should be noted that the estimate of j e comes from the 
fiducial radial profile of the model used to fit the critical 
curve. For small values of 7 e , for example, a critical curve 
can be approximated as: 



(?) = 9, 



0"(fe) 



2 (0"-l + 2 7 e cos(2i/)) (22) 



This naturally means that there is a degeneracy in shape 
such that for fixed shape of a critical curve: 



7 e cxl-0"(0 B ) 



(23) 
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Henceforth, it will be assumed that 7 e is the shear esti- 
mated by fitting the critical curve to a Singular Isother- 
mal Sphere plus an external shear. If we wish to model 
the curve using a different radial profile, a correction of 
7 e (I — tp") must be included. 

We are now ready to Taylor expand equation (|II[) into 
a form that can be estimated only using direct observ- 
ables, and an assumed power-law radial profile for a cir- 
cular lens. In the flux anomaly parameter, only ipn has 
any dependence on the external shear terms. For the 
rest, we may easily relate all terms with a combination 
of radial derivatives and trigonometric functions in rj. 
For example: 



4>i 



ip'(V) sin 2 ^) + e ^"{V) cos 2 (r?) + #o7 e cos(2^) 

To 

and similarly for higher derivatives. Combining all terms, 
and dropping all those quadratic or higher in 7 e , we get: 



.4 



fold ' 



7 e cos(2i/) fl-ip"(6 E ) 



I - ip'' 



n fi-jj" + ie ip'" 

l - V" 



(25) 



Equation (|25j) may seem complicated, but most terms 
are immediately measurable either directly (e.g. 9q) or 
through a fiducial model (77). The shear terms (7 e , v), 
in particular, can be deduced nearly uniquely from the 
observed galaxy and quasar image positions. The various 
radial derivatives can be computed from various circular 
potential models, but assuming that 9q ~ Oe, the range 
of possible flux anomalies is quite narrow. 

Some examples: 

1. Isothermal Sphere (at the Einstein radius) : t/>" = 0, 
ip"' = 0, so: 

.(SIS) _ 7 e cos(2;/) _ rj_ 
f° ld ~ Wo?] 6 

2. Point Source (at the Einstein radius): ip" ~ —I. 
ip'" ~ Aswl/6: 

APS) f cos{2v) 77 
fold ~ ^To 

Note that these are differences of less than two in the 
anomaly parameter, Af Q id, over a fairly wide range of 
potentials. In the next section, we will show that even if 
we use an incorrect global model for a fold lens, we are 
still able to reproduce accurate fold ratios for the smooth 
component of the system. Finally, it will be noted that 
at least one term each of the anomaly parameter models 
has r\ in the denominator. Further, we have noted that 
77 is necessarily a small angle. It is precisely because of 
this form that the anomaly parameter can quite large. 
However, since the external shear appears in the numer- 
ator, and that term is also typically small, we do not 
find any systems in which the estimated fold relation is 
divergently large. 

3. SIMULATIONS 
3.1. Smooth Model Reconstructions 



As a test of the geometric fold approach, we run a num- 
ber of simple model galaxies through lensmodel. The 
source galaxy was chosen to be in near proximity to the 
lens caustic, producing a fold. While most of the the 
source positions were put in "by hand" they were se- 
lected to produce foreground image positions consistent 
with a strict definition of a fold. For us this means that 
the pair needed to be separated by less than half the 
characteristic radius of the system. As KGP point out, 
however, the position along the caustic significantly af- 
fects the expected fold relation. Thus, for our first set 
of simulations, those for a Singular Isothermal Sphere 
with external shear, we've been careful to explore the 
fold relation for those source along the caustic compared 
to those along a radius in the source plane. 

In each case, we have generated a simulated image set, 
including image positions and fluxes, and have assumed 
a knowledge of image parities. We have also assumed 
that the lensing galaxy centroid position is known. For 
the observed images, we then perform a simple lens 
model fitting, assuming only image positions, and using 
a very simple model of a singular isothermal sphere with 
external shear (regardless of the "true" lensing galaxy) . 
We then compare the resulting flux anomaly parameter 
(A f id) to the "true" parameter found from the image 
fluxes. 



3.1.1. Singular Isothermal Sphere + External Shear 

As a first test, the lens galaxy consisted of a simple 
Singular Isothermal Sphere (SIS) with Einstein radius, 
Be = 1, with an additional external shear, 7 e = 0.0.15. 
Though this is a some what higher external shea r than 
typically observed fe.g. lHolder fc Schechterl[2~003| ) . we've 
selected such highly elliptical models as an upper bound 
on reasonable physical systems. Since for all models we 
use an SIS4-external shear model to reconstruct an esti- 
mate of the critical curve, to some degree, this simulation 
is simply a test of our algebra. 

In Fig. ([3]), we plot the measured flux anomaly 
parameters against those estimated by our geometric 
approach. We took a set of sources lying along the edge 
of the caustic, as well as another placed along a radius. 
By far most of the variation in anomaly parameter was 
due to position along the caustic rather than in radius. 
The mean error in the reconstructed anomaly parameter 
is {5 A = A est — A trU e) < 0.02, corresponding to an error 
in the flux anomaly of \SR\ < 0.01. To reduce clutter, 
we have shortened, Af a id, e st to simply A est , here and 
elsewhere. This is far smaller than the factor of ~ 2 
arising from an uncertainty in the radial profile of the 
circular component of the lens. 

3.1.2. Point Source+ External Shear 

As a second test, our true lens consists of a point source 
(again with Be — 1) with an external shear of 0.2 or 
0.3. This is a test that the true radial profile has little 
effect on the reconstructed shape of the critical curve. We 
plot the reconstructed flux anomaly parameters against 
the true anomaly parameter in Fig. 3J Unlike in the 
previous test, points were selected randomly to lie near 
the critical curve such that the observed images satisfy 
the fold condition. 
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Note that (SA) ~ 0.03 for sources near the middle 
of folds, while (SA) ~ 0.04 for sources nearer to cusps, 
resulting in a typical SR ~ 0.02. Since truly anomalous 
flux ratios tend to be in the neighborhood of R ~ 0.5, 
this is the difference of only a few per cent. 



3.1.3. Singular Isothermal Ellipsoid 

In much of our derivation, and in particular, in equa- 
tion (|16[) , we made the assumption that third derivatives 
of the potential contained contributions from only the 
circularly symmetric part of the potential. This is ex- 
actly true for external shears, but not, i n general, for 
ellipti cal potentials. However, as shown bv lKeeton et al.l 
(1997), there is an approximate degeneracy between ex- 
ternal shears and ellipticities, and we will exploit this 
here. As a check of whether our analytic derivation was 
good enough, we have estimated flux anomalies for sim- 
ulated Singular Isothermal Ellipsoids, in this case, with 
a lens ellipticity of 0.5. As above, we've simulated fold 
systems and reconstructed them using the geometric ap- 
proach. Our results are plotted in Fig. [S] 

While cuspy folds are modeled extremely well, with 
(SA) ~ —0.02, more typical folds exhibit a larger sys- 
tematic error, with (SA) ~ —0.09. Of course, this results 
in a systematic error in the measured flux anomaly of 
\5R\ ~ 0.04 for typical fold systems, an effect smaller 
than even typical observational uncertainties for many 
lenses. Since an SIE does not have as many similar 
properties to other mass models that include external 
shear, testing our geometric method with a "true" ellip- 
tical lens re-fit to a simpler mass model with external 
shear should result with values that give some error to 
the flux anomaly. Our results prove that even big differ- 
ences in the properties of the true lens don't contribute 
a significant amount of error to the flux anomaly when 
the assumption is made that the lens is a simple mass 
model with some amount of external shear. 

3.2. The Critical Curve 

In the simulations above, we have shown that the geo- 
metric model can be used to predict the flux anomaly of 
a smooth lens with good precision, regardless of whether 
we use the "correct" model to estimate the critical curve 
of the lens. Because the critical curve must essentially 
thread the observed images, we have contended that the 
shape of the critical curve (and hence the values of i], 7 e 
and v used in equation I25p will largely be independent 
of the smooth lens model used to produce it. 

Starting with an SIS with external shear of 0.2, we 
generated a fold image pair from placing a source near 
the middle of the fold of the caustic. The image positions 
were re-fit using a number of models: 1) A SIS with 
external shear, 2) A point source with external shear and 
3) A Singular Isothermal Ellipsoids. Fig. © illustrate 
the critical curves for the various reconstructions. 

By the shape of the critical curves near the images, we 
have proven that the values of 77, j e and v are indepen- 
dent of the mass model used, thus the values generated 
from the three mass models are generally the same. For 
the models with external shear: (Sv = vsie — vps) — 
0.002 and (Srj = rj SIE - Vps) - 0.020. What is most no- 
ticeable from this plot and the above results is that the 



radii of curvature of all of the external shear models (1 
and 2 in the list above) have very nearly the same radius 
of curvature and orientation. The ellipticity model (3) 
seems to produce relatively different curves and it might 
be expected that they would produce significantly dif- 
ferent models for the flux anomaly parameter. As we've 
seen above in H3.ll this is not the case. This further 
demonstrates that refitting image positions to a very sim- 
ple mass model, such as an Isothermal Sphere with ex- 
ternal shear, still provides a good estimate of the flux 
anomaly parameter even if the true mass model's value 
of 77 based on it's critical curve is not close to the modeled 
values. 

3.3. A Foray into substructure 



The effect of 
quasars has been 
(iMao k Schneider! 



2001 



Dalai k K ochanek 



Williams ct al 



substructure in multiple-imaged 
well-stu died in both sim ulation 
iMetcalf k Madaul 
iDobler k Keetod 
and in a number of o b served 
200 2i: IChiba et al.l 12001 
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20081) ana 111 
systems (IBradac et al.l _ _ 

iMiranda k Jetzeri l2007t iMore et al.l 12009ft . One of 
the major motivations for our geometric approach is 
that we anticipate being able to unambiguously estimate 
substructure within the cluster potential. The shape of 
the critical curve is largely dominated by the smooth 
component of the potential, while the flux anomaly 
depends explicitly on higher derivatives and thus can 
be significantly affected by substructure. Moreover, 
we propose that model-dependent simulations are not 
necessary. 

As a simple proof of concept, we simulated a point 
mass lens with 0e = I and an external shear of 7 e = 0.2. 
We then considered the effect on a fold lens if we placed 
a second point mass lens with 8e = 0.2 (a realistic 4% 
mass perturbation) in the proximity of the fold. In each 
case, the point mass was placed a distance comparable to 
di, and subsequently the substructure mass was rotated 
around the fold pair. Even in an extreme case, the shape 
of the critical curve near the fold - and thus the positions 
of the images - is largely unaffected. Indeed, in almost 
every model (save one), the estimated flux anomaly pa- 
rameter (generated geometrically) was within 20% of the 
anomaly parameter for the unsmoothed system. 

However, the observed flux anomaly was another mat- 
ter entirely. In Fig. [7]we show the dependence of A on the 
angular position of the substructure. Moving forward, it 
is anticipated that substructure can be included in our 
geometric model up to a degeneracy in mass, separation, 
and position angle of the source. 

4. OBSERVATIONAL SAMPLE 

We will now apply our new approach to observed fold 
lenses. We have collected eleven currently known fold 
lenses observed in the radio and IR. A summary of the 
observational references of the systems can be found in 
Table Q] A lens is included as a "fold" if it was identified 
as such in its primary observation paper. In addition, 
wc consider a fold pair if the separation between the two 
is less than Q.59e (where 0e is the best-fit circular lens 
profile) and the next closest pair has separation greater 
than 8e- 

In all the fold lenses analyzed here, we assign image 
"A" to be outside the tangential critical curve and to 
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have positive parity, regardless of the designation origi- 
nally given by the investigators. Likewise, image "B" is 
the image inside the critical curve with negative parity. 

When information about the position of the lensing 
galaxy in a system is unknown (B1555+375, B1608+656, 
B1933+503), we use Keeton's lensmodel software to 
make an estimate of its location. Otherwise, our models 
are based on the positions of the observed images only. 
As discussed in the simulations section (above), the local 
shape of the critical curve (the only piece of informa- 
tion required in our estimates) is not strongly dependent 
upon choice of model. For consistency, we choose to use a 
Singular Isothermal Sphere (hereafter SIS) with external 
shear for the lens. In the case of B0128+437, ellipticity 
is also included to produce a sufficiently accurate recre- 
ation of the observed image positions. 

5. RESULTS 

In Fig. HI we compare our estimates of flux ratio 
anomalies with those found using Monte Carlo analysis 
by KGP, and those observed. A detailed description 
of each system (including three not analyzed by KGP) 
follows. It is noteworthy that most of our reconstruc- 
tions produces similar estimates to flux anomalies as 
those found by KGP. However, many systems exhibit 
fold relations inconsistent with smooth lensing models. 
For lobe-dominated radio sources, this likely suggests 
substructure within the lens. For core-dominated 
sources, microlensing may be at work, depending on 
the beaming of the radio lobe. In the latter case, this 
ambiguity could presumably be resolved by exploring 
the lenses in the time domain. 

B0128+437 is a system with both an extremely high 
flux anomaly, as seen at several diff erent frequencies 
([Biggs et al.ll2004t iPhillips et aT1l2000h . as well as mul- 
tiple components in the source. Images A, C and D have 
been resolved into three different components embedded 
into a more extended jet. Since image B's components 
are not well-defined, we c an only directly consider the 
integrated flux ratio here (jBiggs et al.ll2004h . Moreover, 
since the observed components of images A and B are 
clearly resolvable, shap e analysis naturally w ould yield 
additional information (|Koopmans et al.ll2003|) . 

Observationally, this system exhibits a very large flux 
anom aly in the range of R = 0.263 (|Koopmans et al.l 
12001 to R = 0.582 (jBiggs et al.l 12004ft . KGP have 
studied this system and found from simulations an 
expected flux anomaly in the range of —0.1 to 0.4, with 
a preferred value around R = 0.25. Our own analysis of 
this system predicts a value of R — 0.161 for an SIS+ES 
and R — 0.176 for a PS+ES, well within the range pre- 
dicted by simulation. Further, high-resolution imaging 
of the system resolves unambiguous lobes, suggesting a 
source too large to be affected by microlensing. This 
points unambiguously to additional substructure in the 
lensing galaxy. We will explore the question further 
constraining substructure in future work. 

B071 2+472 is classified by KGP as both a cusp and a 
fold. B0712+472 violates the cusp relation, as discov- 
ered in Keeton et al. (2003), but only in optical data, 
which is evidence more for microlensing than for small- 
scale structure. The four images are relatively circular 



and easily distinguishable from each other, but there is 
a faint "bridge-like" feature betwe en components A an d 
B, which are very close together ( Jac kson et al.l [T998f) . 
suggesting a slightly extended source at the highest reso- 
lution. Spectroscopically, the quasar has a relatively flat 
spectrum which, along with most photometric estimates, 
imply a core-dominated source. 

The KGP Monte Carlo simulations predict Rfdd to be 
between -0.1 and 0.15. Our value for Rfoid from fitting 
the lens to an SIS+ES is -0.023, in agreement with 
the KGP simulation s . The observed radio observ ations 
(jJackson et al.l 119981 120001: IKoopmans et al.l I2003D also 
agree with the predictions. 

B1555+375 images A and B appear in radio observa- 
tions as one extended object with two separate peaks of 
brightness. The fourth image was not initially observed 
but was predicted, searche d for and consequently discov- 
ered (jMarlow et al.lll999bh . KGP find an expected range 
for the value of Rfoid for this lens to be from about 
-0.07 to 0.2. Our model using an SIS+ES produces an 
Rfoid of about 0.003, which goes up to 0.062 if we use a 
PS+ES. The observed Rfoid, which is greater than 0.23 
in all radio observations, exceeds our estimate, suggest- 
ing the system has significant small-scale structure or 
potentially microlensing. The imaging has insufficient 
resolution in the radio to determine whether the source 
is extended. However, no double-lobes or jets are clearly 
visible. 

B1608+656 has a core-dominated source with a com- 
plex lensing struc ture in the form of a second galaxy. 
IMvers et all ()1995f ) claim that B1608+656 is simple com- 
pared to other lens systems because they were able to re- 
produce the image configuration using a single, elliptical 
lens galaxy in their model. The ellipticity of the lens in 
their models turns out to have been a sufficiently accu- 
rate recreation of the combination of lensing galaxies in 
the real system. 

KGP predicts Rfoid to be in an unusually small 
range, from 0.3 to 0.43. Observed Rfoid values, 

which range from 0.327 to 0.51 6, m atch up well 

(IMvers et aLl 119951: IFassnacht et all 119991: iSnellen et all 
119951: iFassnacht et al.ll2002fl . Our estimate for a SIS+ES 
is R fold = 0.361, and for PS+ES is R fotd = 0.577, 
which is higher than both observations and other pre- 
dictions, suggesting a more isothermal mass distribution. 

B1933+503 has ten images as a result of three lensed 
images. Two of the sources result in quad formations 
and the third appear s doubly imaged. As noted by 
IKochanek et al.l (|2004D . the three sources seem to cor- 
respond to a core and both radio lobes. The four images 
that appear in a standard quad formation are a cross-like 
fold image set. 

KGP model Rfoid in a rather wide range from about 
-0.1 to 0.45, with a sharp peak at 0.05 and a shorter, 
wider peak at 0.3. The observed flux ratios (at any 
wavele ngth) are much hi gher: in the range of 0.580 to 
0.722 (|Svkes et al.l 119981 ). Our geometric estimate of 
the flux ratio using a standard SIS+ES is —0.233. This 
anomalous outlier presents a limit to our approach, since 
iNairl (|1998ft model with a lens ellipticity of 0.81. 
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B1938+666 is a lens with two sources, quite possibly ex- 
tended radio lobes. The radio lobe-dominated emission 
picture is further supported by a spectral index of about 
0.5. One source creates a four image fold configuration, 
which we use for further analyses. The other source pro- 
duces a double (jKing et al.lll997l) . In the near- infrared 
and optical wave bands, a slight Einstein ring appears 
(IKing et al.lll998h . 

The two fold images have a very small flux anomaly 
as detected using ME RLIN, but a muc h higher value, 
-0.436, using VLBI ((King et alJll997t ). A geometric 
reconstruction of this lens is much closer to the high 
value, which either suggests substructure or perhaps sig- 
nificant time variability. The observed images are re- 
solved enough to measure shapes, which would poten- 
tially further constrain mass models. 
HS0810+2554 has not been observed in the radio, 
but we used the near-infrared CASTLES observations 
which yield a flux ratio anomaly of 0.2 74. The lens 
was discovered by iReimers et al.l (|2002f ). who noted 
that this system's configuration is very similar to 
that of PG1115+080. Our geometric analysis yields a 
value much closer to zero, however, suggesting either 
small-scale structure, or possibly differential absorption. 

MG0414 + 0534 has been observed over a wide range 
of frequencies and epochs. From optical observations, 
there is a non-distinct arc that passes through the fold 
ima ge pair, Al and A2, an d the next closest image, 
B (En gonin-Willaime et alJ ll999). Radio data clearly 
shows the resoluti on of the images in th is fold configura- 
tion lens system (|Katz fc Hewittl[T993fl . 

MG0414+0534 has lobe-dominated emission and an 
extremely steep radio spectrum. It also has a satellite 
galaxy, "object X", that is not particularly close to the 
fold image pair. In order to match the image positions, 
KGP needed to include object X in the model. The 
Rfoid values estimated from KGP range from about 
0.0 to 0.2 with a peak at around 0.07. Observed data 
has Rf n u values r anging from about . 047 to 0.073 
(IKatz et all 119971 : IKatz fe Hewitt] 119931 : Hewitt et all 
11993 and CLASS) which agrees our estimates. We did 
not include object X in our analysis, which resulted in 
a somewhat lower flux anomaly ration than observed or 
estimated by previous analysis. 

MG2016+112 has two separate sources, lensi ng to a quad 
and a double, to produce a total of six images ([King et al.l 
1997). Not all ima ge components are vis ible in the vari- 
ous wavelengths. iSchneider et al.l (|1985l) observed three 
images and one lens galaxy in both radio and opti cal ob- 
servations. One year later, ISchneider et al.1 (fl986) found 
an additional image. CASTLES has observed three im- 
ages and four lens galaxies. 

We reconstructed this lens with a lens position from 
the CASTLES op tical data, and image positions from 
IMore et al.l pool 's 1.7GHz data. They conclude, 
"there is no significant substructure or any other effects 
that might affect the flux densities of the images." Our 
reconstruction yields a flux fold ratio consisten t to that 
observed, further supporting IMore et al.l (|2009l ) . 

PG1115+080 has thus far been measured only in the 
mid-IR. The NICMOS images of PG1115+080 suggest a 



quasar host galaxy that lenses to an Einstein ring around 
the four fold configuration images. The lensing galaxy 
has no substruct ure, and there is no blatant flux anomaly 
in the infrared (|Impev et al.lll998h . 

KGP predicts the value of Rfoid to range from -0.05 to 
0.25, with a peak at about 0.1. This is consistent with 
observed values (|Chiba et al.l 120051: llmpev et all 119981: 
IChiba et al.ll2005h . Our own reconstructed fold relation 
is comfortably in this range as well. 

From our reconstruction of the image positions using 
an SIS+ES yield an Rfoid value of about 0.045 which 
fits our observed infrared data from IChiba et all (l2005h 
and fits perfectly in the range determined by KGP. This 
further supports the thought that this lens has little to 
no substructure. 

SDSS1004+4112 is a very massive system (jlnada et al.1 
l2003h which has only thus far been reliably measured in 
the near-IR. It contains a five-image fold-configuration 
system with a complicated galaxy cluster at the heart of 
its lens (|Inada et al.ll2005D . 

KGP predicts SDSS1004+4112 to have R fold values 
between and 0.25, with a sharp peak at 0, with 
higher flux ratios less likely. Observations fall toward 
the high end of this range, as does our own geometric 
reconstruction. This system is known to have complex 
lens structure, but no flux ratio anomaly occurs in the 
fold pair images. 

SDSS1330+1810 is an apparently lobe-dominated source 
observed in the near-infrared and optical wavebands with 
the Magellan, UH88, and ARC3.5m telescopes. The 
fold relation is vastly different in the optical and near- 
infrared, thus lOguri et al.l (|2008l ) predicts dust to be 
the cause of this anomaly. From analyzing the images, 
there may be some substructure in t he form of a cluste r 
of galaxies positioned near the lens (|Qguri et al.ll2008D . 

We performed a reconstruction of this lens with only 
the near infrared data in the J, H 7 and K s bands. The 
values of Rfoid from observed near-infrared data ranged 
from 0.101 to 0.151. Our geometric reconstruction pre- 
dicts a smaller fold relations, of only 0.007, suggesting 
contribution of substructure. 

6. CONCLUSIONS AND FUTURE WORK 

Thus far, we have focused on estimating the flux 
anomaly for fold lenses using entirely geometric argu- 
ments. This is approach is useful because it produces 
robust estimates of flux fold relations without the need 
to use Monte Carlo simulations. In fact, the only model- 
specific, non-observed, parameter to be adjusted is the 
index of the radial profile which can be tuned analyti- 
cally. This is both easier to apply than past approaches, 
and also provides a much more intuitive understanding 
of the underlying structure of the lens. However, it is 
worth noting that Monte Carlo analysis is still necessary 
to get a deeper understanding of the underlying variance 
in the fold relation distribution. 

This work is an important first step in further geo- 
metric, semi-analytic analysis methods, and has a num- 
ber of potenti a l futu re offshoots. For one, following 
ICongdon et all ((2008), we will extend our analysis to in- 
clude analytic estimates of time delays. Since time delays 
are functions of only the potential itself, the analysis is 



expe cted to be significantly more robust. Further, in 
we performed a tentative analysis of a system with 
substructure. Future work will be required to identify 
the degenerate families of substructure which can identi- 
fied by geometric analysis. Another logical next step is to 
perfo rm a geometric analysis of cusp lenses (jKeeton et al.l 
2003). The central focus of this work, however, is based 
on the premise of uniquely identifying substructure in 
lenses and our long term focus will be to provide a sys- 
tematic analysis of substructure in quad lenses using this 
framework. 
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visible. Errorbars are listed for all observations or papers where errors are 
reported. All distances and Einstein radii are in arcseconds. Flux ratios are 
dimcnsionlcss. Estimates of the Einstein radius are taken from KGP. The 
geometrically modeled flux ratio anomalies are fit using a Singular Isothermal 
Sphere plus external shear profile. 
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Fig. 1. — A typical fold lens. Left: Source plane with caustic curves. Right: Image plane with critical curves. Images with negative 
magnification (parity) are denoted with filled circles, while images with positive magnification (parity) are open circles. Throughout the 
paper (and in most papers on fold lenses) Image A is the positive parity fold image, and Image B is the negative parity fold image. Following 
KGP, di is the distance separating A and B. 
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Fig. 2. — The relative orientation of the fold images (A,B), the lens center, and the critical curve in our model. The point, V, represents 
the position on the critical curve where the two images would meet if the source were moved closer to the caustic in the source plane. It is 
also the point around which all derivatives of the potential are defined. 
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Fig. 3. — The measured and estimated flux anomaly parameter for a Singular Isothermal Sphere (SIS) lens with an external shear. This 
simulation, and those in subsequent figures, was produced using the lensmodel package. In each case, the lens has an Einstein radius of 
8 e = 1, and in this model, there is an external shear of 7 e = 0.15. As KGP have shown, we get dramatically different flux anomalies 
depending on where the source image lies along the caustic. Triangles represent sources near the edge of a caustic. We have uniformly 
placed 18 sources approximately 10% the characteristic scale from, and along the edge. The large flux ratios arise for sources nearer to the 
cusps. Unsurprisingly, those sources near cusps more strongly resemble "cusp lenses." Since cusp lenses are expected to have a different flux 
relationship (in which the bright image equals the sum of two dimmer ones), it is unsurprising that the flux anomaly between two images 
can be rather large. Indeed, this is precisely the result found by KGP. In every case, the geometric reconstruction technique produces a 
very good fit ((SA) = 0.02, SR < 0.01). The 8 diamonds represent images placed radially inward from the center of the caustic to within 
10% of the center of the source. Those with a small (SA ~ 0.1) deviation from prediction are those closest to the center, and most closely 
resemble quad lenses, rather than true folds. The diagonal line is meant as a guide and has a slope of 1. 
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Fig. 4. — The measured and estimated flux anomaly parameter for a Point Source (PS) lens with an external shear. In each case, the lens 
has an Einstein radius of 6 e = 1, and various systems were modeled with an external shear of 7 e = 0.2 and 0.3. Since most of the variation 
in flux anomaly arises for variations along the edge of the caustic, we use a slightly different set of symbols in this figure and subsequently 
than in the previous one. Triangles represent sources near the edge of a caustic, while squares represent sources closer to cusps. In both 
cases, there is a small systematic error in the anomaly parameter of (8 A) ~ 0.08. This corresponds to a typical error in the flux anomaly 
of 8R ~ 0.04. 
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Fig. 5. — The measured and estimated flux anomaly parameter for a Singular Isothermal Ellipsoid (SIE) with an ellipticity of 0.5. The 
fiducial lens has an Einstein radius of 8 e = 1. As in the previous figures, triangles represent sources near the edge of a caustic, while squares 
represent sources closer to cusps. Cuspy folds seem well modeled, with (SA) = —0.02. However, sources further from the cusps had a larger 
systematic error, with (SA) ~ —0.09. This one-sided offset is primarily due to a flaw in our assumptions. Throughout, we've assumed that 
the shear only contributes in second derivatives of the potential. In future work, we may relax this condition somewhat and allow for an 
intrinsic elliptical potential. 
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Fig. 6. — The fold lens images for a simulated quad system with a Singular Isothermal Sphere (Bp = 1) and an external shear with 
7e = 0.2. Only the two fold images are shown in this plot, but the critical curves are fit to all four image positions, and the image centroid. 
Curves are fit for several models with external shear, as well as potentials with intrinsic cllipticity. The solid line critical curve represents a 
Singlular Isothermal Sphere with external shear, the dashed line represents a Singular Isothermal Ellipsoid and the dotted line represents 
a Point Source with external shear. 
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Fig. 7. — As described in the text, the observed flux anomaly parameter for a variety of simulated fold lens systems. For each, a point 
mass with external shear was placed at the origin. In addition, a substructure with 6e = 0.2 was placed a fixed separation (approximately 
equal to the image pair) such that the angle between the image separation and the secondary mass was S. 
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Fig. 8. — A comparison of our simulated Isothermal Sphere+external shear estimates of the fold relation (x-axis), with those estimated 
by KGP (diamonds with errorbars) from Monte Carlo simulation, and observed (triangles). This sample consists of the 10 quad systems 
in the radio and IR for which analysis was done by both groups. Note that our results have a very close correspondence ((AR) ~ 0.02) 
with those done using Monte Carlo simulations, but make much simpler assumptions. Further, for about half the systems, the observed 
flux anomaly is very similar to those estimated by both groups. 



